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ABSTRACT 


In 2009, a monthlong study was conducted in the Western Philippine Sea. This 
experiment collected oceanographic and acoustic data from two moorings, the 
Distributed Vertical Line Array (DVLA) and T1 moorings. The T1 mooring used a 
1,095m depth source to send broadband acoustic signals to DVLA receivers 
located between the depths of 800m and 1,385m. Jacob A. Fischer conducted 
the internal tide analysis presented in Chapter II. The oceanographic data were 
used to quantify the characteristics of internal tides in the sample area. Spectral 
analysis and mode fitting, shows that most of the energy within the internal tides 
is due to mode 1 of the tidal frequencies. Weston R.T. Coby conducted the 
acoustic analysis presented in Chapter III. The acoustic data were first used to 
identify acoustic timefronts across the observed depths, and track the timefronts 
arrivals. A ray prediction model was then used to match eigenray properties to 
the observed timefronts. Finally the Viterbi program was used to track timefront 
arrivals, and intensity for statistical analysis. Weston R.T. Coby and Jacob A. 
Fischer jointly compared the internal tide and acoustic variability presented in 
Chapter IV, proving internal tides significantly impact sound speed fluctuations in 
the Philippine Sea. 
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I. INTRODUCTION 


Weston R.T. Coby, and Jacob A. Fischer jointly prepared the introduction 
chapter. 

A. NAVAL IMPORTANCE 

To maintain the upper hand in the future subsea battle space, the Navy 
must gain a better understanding of how oceanographic processes affect the 
ability to detect and identify potential threats. Currently, acoustic detection and 
resolution methods are favored in the subsea environment. Sound speed 
perturbations in the ocean medium adversely effect object classification and 
resolution capabilities. Internal tides are a major source of sound speed 
perturbations in the open ocean. Identifying and classifying the relationship 
between sound speed perturbations due to internal tides and acoustic travel 
times and intensities will allow for future advancements in ocean acoustic 
technology. Understanding the ocean acoustic variability improves our target 
identification and detection capabilities. 

B. PHILIPPINE SEA 

The Philippine Sea is an oceanographically significant marginal sea in the 
Pacific Ocean, northeast of the Philippine islands. It is considered one of the 
most dynamic oceanographic environments in the world, and it is an area that 
generates very strong internal waves and tides. Several faults and fracture zones 
border the Philippine Sea, and its basin is host to island chains (Philippine, 
Marianas, Ryukyu), trenches (Philippine, Mariana) and ridges making it a region 
of geographical significance (Figure 1). 
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Figure 1. Philippine Sea Topographic map showing location of significant 
features: trenches, archipelago, island chains. From 
http://mapsof.net/map/philippine-sea-location#.Uf89EVOf9IU 

Internal tides are internal waves propagating at tidal frequencies. Internal 
tides (baroclinic tides) are generated by the barotropic tidal forcing of water over 
significant bathymetric features. The more significant the bathymetry and 
barotropic tidal forcing within the region, the more intense the internal tides 
(Garrett and Kunze 2007). The Philippine Sea is host to several potential 
generation sites for internal tides due to its complex and rough sea floor. The 
Luzon Straight in particular is an area known to generate strong internal tides, 
but their generation at the site is variable (Colosi et al. 2012). Several 
experiments have been conducted to study the characteristics, and propagation 
of internal tides in the northern Atlantic, central Pacific, and eastern Pacific, yet 
few have been conducted in the western Pacific. Figure 2 shows the modeled 
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displacement (flux) - dominated by mode one- of the internal tides. The 
significance of the Philippine Sea is best seen in Figure 2 by the relatively large 
displacements in the northern Philippine Sea as compared to the north central 
and northeastern Pacific Ocean. 



-20 -15 -10 -5 0 5 10 15 20 

Normalized interface displacement (m) 

Figure 2. Shows the propagation of internal tides as measured through 
displacements (After Garrett and Kunze 2007) 

C. 2009 NPAL PHILIPPINE SEA EXPERIMENT 

The North Pacific Acoustic Laboratory (NPAL) is a term that was first used 
to describe a series of acoustic propagation experiments using low frequency 
broadband sources, and various large aperture vertical and horizontal receiving 
arrays in the mild (low-level eddy variability) northeast and north central Pacific 
Ocean from 1989-2005. From these initial experiments, the background sound- 
speed field, small-scale sound-speed variability due to internal waves, spice, and 
noise sources were determined (Worcester and Spindel 2005). A second series 
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of deep-water acoustic experiments were conducted from 2009-2011 in the 
oceanographically and geologically complex region of the Philippine Sea. 



Figure 3. NPAL Philippine Sea Experiment area (From Worcester 2008) 


The series of experiments consisted of the NPAL pilot study/engineering 
test (Spring 2009), NPAL Philippine Sea Experiment (2010-2011), and Ocean 
Bottom Seismometer Augmentation of the NPAL Philippine Sea Experiment 
(2010-2011). The Philippine Sea experiments utilized moored and ship 
suspended low frequency transmitters to send signals to the newly developed 
Distributed Vertical Line Array (DVLA), acoustic Seagliders, bottom 
seismometers, and a towed hydrophone array to record the transmissions and 
ambient noise. Some of the goals of these experiments were to better 
understand (Worcester 2013): 

• The impacts of ocean processes on acoustic propagation. 

• The ability to generate a useful ocean model that can aid in 
improving acoustic predictions. 
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• Scattering, as it relates to spice and internal waves. 

• How ambient noise varies in time and space. 

• The connection between the seismic field in the seabed and sound 
propagation in the water column (Worcester 2013). 

This thesis focuses on the data collected during the 2009 NPAL pilot 
study, conducted between April-May 2009, and sponsored by the Office of Naval 
Research (ONR). The monthlong pilot study was conducted using a single 
transmitter to test the DVLA in preparation for the follow-on yearlong experiment, 
which used six transmitters. The pilot study provided initial acoustic observations 
in the northern Philippine Sea, and verified the proper function of the equipment 
used in the experiment the following year (Worcester 2008). 

During the pilot study acoustic and oceanographic data were collected 

from two moorings (DVLA and T1) within the region, and from multiple 

conductivity-temperature-depth (CTD) profiles. Broadband signals were 

transmitted from the T1 mooring (Number 1 in Figure 3), travelling a horizontal 

distance of 185 km, to the DVLA mooring (as shown in Figure 3). The NPAL pilot 

study is unique because this is the first time that scientists have been able to 

track timefronts throughout the water column in the Philippine Sea. Through the 

NPAL pilot study and its use of the DVLA, the scientific community is gaining the 

ability to track and analyze acoustic timefronts throughout the ocean column, in a 

complex environment to better determine the impact of environmental features 

such as internal tides on acoustic propagation. The acoustic information gathered 

required several layers of filtering in order to quantify the physical properties and 

statistics for each acoustic ray at specific depths. The signal-to-noise ratio (SNR) 

was calculated from the acoustic data collected in the pilot study, and used to 

isolate and track the acoustic ray timefronts and the travel time fluctuations over 

the duration of the experiment. Figure 4 is a plot of the observed timefronts for a 

single transmission collected from the DVLA receivers. Figure 4 shows the 

pattern of a wavefront sweeping by the DVLA receivers at a fixed range, 

displaying the timefront arrival times, and timefront variation with depth. 

Calculating SNR enabled us to separate the transmission acoustic timefronts 
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from the background noise (blue background). Tracking the individual timefronts 
across depth is important because it allows us observe and analyze the evolution 
of sound rays. The arrival times for the individual timefronts vary significantly with 
increasing depth as shown by Figure 4. 



Figure 4. Timefront plot displays the separate eigenray arrival times for DVLA 
receivers between 800m and 1,385m for a single transmission. 


Figure 5 is a plot of the timefront receptions measured at the 1,000m 
DVLA receiver over the duration of the experiment (y-axis). Since Figure 5 only 
uses the SNR data from the receiver at 1,000m it shows how individual 
timefronts vary in time at a fixed depth. 
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1,000m Arrivals 


40 



Travel time (seconds) 


Figure 5. Time series plot of timefront arrivals at the 1,000m DVLA receiver 
shows arrival time fluctuations for individual timefronts throughout 
duration of the pilot study. 

Moving on to the oceanographic data, there were temperature, salinity, 
and pressure measurements taken at both moorings. These data were exploited 
to track isopycnals (surfaces of constant potential density) and isotherms 
(surfaces of constant temperature) between 120 and 2,000 meters (Figure 6), 
and to calculate sound speed fluctuations stemming from internal tides and 
intrusive oceanic structures. 
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Internal Wave Displacements at DVLA Internal Wave Displacements at T1 



04/16 04/17 04/16 04/17 


Date (2009) Date (2009) 

Figure 6. Time series of internal wave displacement by depth at DVLA and T1. 

Derived from tracked isopycnals and isotherms (After Colosi et al. 

2012 ). 

Oceanographic variability (Figure 6) and the acoustic variability (Figures 4 
and 5) are related to each other by the sound speed changes caused by ocean 
processes. Of these processes, internal tides and random internal waves 
contribute the most to sound speed variability. Vertical tidal displacements within 
the water column cause sound speed variability due to their interaction with 
acoustic signals propagating through the ocean. As the acoustic signal 
propagates, parcels of water within the tidal displacements interact with the 
acoustic signals. These parcels of water are physically different than the parcels 
of water that they have displaced. The differing physical properties (temperature, 
salinity, and pressure) of the shifting parcels of water will either slow or speed up 
the acoustic signals passing through them. The potential sound speed gradient is 
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used in determining sound speed perturbations because it is a good marker of 
water column stability as it relates to internal tides (Munk et al. 2009). Equation 
1.1 is used to relate root mean square (RMS) tidal displacement to sound speed 
perturbation (Colosi et al. 2012). 


8 c = £ 


dc 

dz) 


8c- Sound speed perturbation 
Vertical displacement (RMS) 

f dc^ 


v dzj 


p - Potential sound speed gradient 


( 1 . 1 ) 


The Garrett-Munk (GM) spectrum describes the random internal gravity 
wave field (Munk 1980). The GM spectrum is significant because it discovered a 
frequency-wavenumber distribution for internal waves that appeared to be 
“universal” in space and time, and is “often used as a representative statistical 
description of the internal wave field (Levine, 2002, p. 3166).” Diurnal and 
semidiurnal internal tides which are not included in the GM model show energy 
levels which are approximately equal to those of random internal waves. The 
random internal wave energy measurements at the T1 and DVLA moorings were 
relatively constant, so the diurnal and semidiurnal internal tides were separated 
by their energy fluxes at tidal frequencies (Colosi et al. 2012). 


D. PURPOSE 

The analysis presented in this thesis serves to address three main areas: 
[1] characteristics of observed internal tides, [2] acoustic variability occurring at 
tidal frequencies, and [3] correspondence between moored observations of 
internal tides and tidal acoustic variability. This thesis seeks to enhance our 
understanding of the subsea environment and to improve our understanding of 
how oceanographic processes influence visibility in the sub-sea battle space. 
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II. INTERNAL TIDE ANALYSIS 


Jacob A. Fischer conducted the internal tide analysis presented in this 
chapter. 

A. LINEAR WAVE THEORY: EQUATIONS OF MOTION AND MODES 

In deriving the equations of motion (Equations 2.1-2.3, 2.6, 2.7), we 
assumed the ocean to be an incompressible fluid and used the Boussinesq 
approximation; assuming that changes in density throughout the water column 
were negligible except as they related to terms associated with the force of 
gravity. (Stewart 2013). 


du 1 dp 

dt p 0 dx 

(2.1) 

dv 1 dp 

— + fu- -— 

dt p 0 dy 

(2.2) 

dw 1 dp p 
dt p 0 dz p 0 

(2.3) 


The first three equations of motion (Equations 2.1-2.3), the momentum 
equations, stem from Newton’s second law of motion, which governs the 
response of a parcel of fluid’s momentum as a result of forces acting on it; those 
forces being pressure p, the Coriolis force /, and gravity g. Thinking about the 
ocean in terms of a Cartesian coordinate system, we get forces acting on a fluid 
parcel and the resulting motion in three-dimensional space. The three- 
dimensional space, with respect to the ocean, is broken down into three 
directions of motion and their velocity counterparts. The east/west direction and 
velocity are defined by * and u , north/south by y and v, and depth by z and w. 
Within the x-y plane, forces governing the momentum of the fluid are 
predominantly pressures acting on the fluid ( p) and the force exerted as a result 
of earth’s rotation, the Coriolis force (/). 
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In order to arrive at Equation 2.6, we needed to apply the continuity 
equation (Equation 2.4) and the incompressibility relation (Equation 2.5). The 

N 2 S_<*£ 

Buoyancy frequency squared is equal to p ■ 


^+ u Vp + pVw =0 
at 


Vm =0 


d_ 

dt 


g- 






= wN z 


(2.4) 

(2.5) 

( 2 . 6 ) 


From the principles behind the conservation of mass, for a parcel of water 
moving within a fluid at velocities significantly less than the speed of sound, the 
density (mass/volume) is assumed to be a constant value. So, if the density 
doesn’t change, i.e. the particle does not expand or contract, the divergence of 
the velocity vector field must be equal to zero. This is the concept behind the 
incompressibility relation (Equation 2.7) (Stewart 2013). 


du + dv + dw _ 
dx dy dz 


(2.7) 


The equations of motion were used to derive the linear internal wave 
equation for vertical displacement E,(r,t) (Equation 2.8). Important information for 
the wave equation is the local Coriolis parameter / and the buoyancy frequency 
profile N(z). Figure 7 shows the N(z) profile obtained from the NPAL pilot study. 


Y a 2 Y 

r d 2 d 2 ' 


f a 2 .Y 

a 2 " 

\ dt ) 

v dx 2 dy 2 j 

+ 

v +/ 

l J 

az 2 


Z(r,t) = 0 


( 2 . 8 ) 
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Figure 7. NPAL pilot study buoyancy frequency profile (After Colosi et al. 2012). 

The modal solution in Equation 2.9 is significant because it has within it, 
the components of vertical (standing wave) and horizontal (travelling wave) 
motion. The vertical motion is governed by the mode functions, \\/ (z) and the 

horizontal motion is governed by the propagation term, e {,{kx+l} wt)) . Combined, 
they describe wave motion as vertical displacement over time, and the 

nature of internal tide vertical scales. 

^(r,0 °c \| f .( z )e {l{kx+ly - m)) (2.9) 

Two examples of the solution to the linear internal wave mode equation 
(Equation 2.10), Modes 1 and 10, at the diurnal and semidiurnal frequencies are 
shown in Figure 8. In Equation 2.10, the horizontal wavenumber eigenvalue k j is 
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equal to ^(k 2 +1 2 ). The diurnal and semidiurnal tidal constituents (frequencies) 
are represented by co in Equations 2.9 and 2.10. The diurnal frequencies 
correspond with roughly one high and low tide per day whereas the semidiurnal 
frequencies correspond to two highs and lows. It can be seen from Figure 8 that 
there is frequency dependence with respect to the wave modes. 


d2 Vj 

dz 2 



f 

V 


N 2 - to 2 ") 

r .2 f 1 

« ~f ) 


V;= 0 


( 2 . 10 ) 


Mode 1 for internal tidal frequencies Mode 10 for internal tidal frequencies 




Figure 8. Internal wave modes 1 and 10 for diurnal and semidiurnal frequencies. 

This shows mode shapes have a dependence on frequency. Though, 
when comparing semidiurnal and diurnal frequencies, the difference is 

minimal. 


By solving the wave equation, we can use the wavenumbers to get the 

horizontal wavelength (Table 1). The wavelengths X in Table 1 are broken down 
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by mode number and tidal constituent (01, K1, PI, etc.). The wavelengths are 
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related to wavenumber k by the following relation, X = —. 

k 



Wavelength (km) 

Diurnal Constituents 

Semidiurnal Constituents 

01 

K1 

PI 

Q1 

M2 

S2 

N2 

K2 

Mode 

# 

1 

516 

433 

438 

576 

167 

161 

171 

160 

2 

246 

206 

209 

274 

79 

76 

81 

76 

10 

49 

41 

42 

55 

16 

15 

16 

15 

20 

24 

20 

21 

27 

8 

8 

8 

8 

50 

10 

8 

8 

11 

3 

3 

3 

3 


Table 1. Wavelengths for diurnal and semidiurnal constituents, and their 

respective mode number. 

As can be seen from Table 1, the mode one horizontal wavelength of 
the principle lunar semidiurnal tidal constituent (M2) is around 170km. This is 
comparable with previous work done by Kerry et al. (2013). Preliminary 
analysis conducted in this section and the theory involved was critical to the 
analysis of internal tide energy (Section ll-C) and in generating a regional 
model (Section ll-E). 

B. DIURNAL AND SEMI-DIURNAL 

The tracked isopycnals and isotherms at each mooring served as an 
accurate representation of displacement. The displacement time series used in 
this thesis were previously filtered by Colosi et al. (2012) to remove frequencies 
greater than 500 cpd. To verify the strong diurnal and semidiurnal signals present 
at both DVLA and T1 during the monthlong experiment, a spectral analysis was 
conducted on the tracked isopycnals for several depths at both locations. 

At each depth, the mean displacement was subtracted from the data to 
find the displacement amplitude about the mean. A Hanning window was applied 
to compensate for leakage due to the time series not being periodic and to filter 
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out some of the ultra low frequencies. The data were converted from the time to 
the frequency domain by taking the discrete Fourier transform. Sampling 
frequency was 1/60 seconds. 

To find the peak frequencies, the amplitude spectra were averaged across 
depth. Results showed that below the surface layer, the internal wave structure 
between DVLA and T1 is dominated by waves traveling at the diurnal and 
semidiurnal frequencies (Figure 9). Since we were concerned with quantifying 
the effects of internal tides on acoustic variability, only the diurnal and 
semidiurnal frequencies were isolated in further processing of the data. 
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Figure 9. Average of amplitude spectra across depth highlights peak amplitudes 
at diurnal and semidiurnal frequencies. The yellow vertical line 
corresponds to the intertial frequency. The blue vertical lines 
correspond with diurnal frequencies (Oi,Ki,Pi,Qi). The green vertical 
lines correspond to semidiurnal frequencies (M 2 ,S 2 ,N 2 ,K 2 ). 
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At DVLA, the diurnal frequencies dominate but at T1, the semidiurnal 
frequencies dominate (Figure 9). This is not an unexpected result due to the 
complex bathymetry at the tidal generation zone contributing to spatial 
interference (Colosi et al 2012). Upon determining that the diurnal and 
semidiurnal frequency bands were the most dominant in the internal wave 
structure, a least-squares fit in time was applied to the displacement data (Figure 
10). A weighted least-squares method (Equations 2.11-2.12) was used in fitting 
the unfiltered data to 8 diurnal frequencies (Oi,Ki,Pi,Qi,Ji,NOi,ALPi,OOi) and 7 
semidiurnal frequencies (M 2 ,S 2 ,N 2 ,K 2 ,L 2 ,MU 2 ,A 2 ) (Stewart 2013). 


P=[E r R 1 E ]“ 1 

A ffl =PE r R“\. C) 


( 2 . 11 ) 

( 2 . 12 ) 


From Equations 2.11 and 2.12, P j s a matrix of tidal amplitude 
uncertainties, E js a matrix that contains the estimated internal tidal frequencies, 
R is the squared data error and Ap,°o is a matrix of tracked internal wave 


amplitudes (isotherms and isopycnals). Weighted amplitudes A m from the least- 
squares fit were multiplied by the sine and cosine of the tidal frequencies in each 
of the 2 frequency bands. Then, the components for each frequency were 
summed up with respect to time in order to create the filtered wave spectrum for 
each frequency band. Time sampling in data processing was every 30 minutes. 
Two separate time series (diurnal and semidiurnal frequencies) were generated 
for each depth. Data error in performing the tidal fit was assumed to be ,6m, as 
determined by Colosi et al (2012). 
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Diurnal Least Squares Fit at DVLA 



Semidiurnal Least Squares Fit at T1 



Figure 10. Least-squares fit to displacement for diurnal tidal components at DVLA 
and semidiurnal tidal frequencies at T1. Raw data is red and blue. Fit 

data is black. 


C. INTERNAL TIDE ENERGY 

Internal tide displacement fluctuations are dependent upon depth, as can 
be shown by observing the RMS displacements for the tidal components at both 
moorings in Figure 11. It is apparent from Figure 11 that the RMS displacement 
does not follow a strict mode one shape across depth, inferring that the 
waveform is more complex than the mode one structure. Applying Equation 1.1 
to the RMS displacements, we were able to get sound speed perturbations 
(Figure 12). 
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Depth (m) 


RMS Displacement at DVLQ 


RMS Displacement at T1 



Figure 11. RMS displacements of tidal components at both moorings. 
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Potential sound speed gradient 


Sound speed perturbation 




Figure 12. Potential sound speed gradient and sound speed perturbations for T1 

and DVLA tidal components. 

To determine the contribution of different modes to the structure of the 
waveform, we had to determine the energy contributions from these modes 
(Colosi et al. 2012). 50 modes were fit to the semidiurnal and diurnal 
displacement time series at both moorings using Equations 2.13 and 2.14, and 
the Gauss-Markov estimator (Equations 2.15-2.16) (Colosi et al. 2012). In 
Equations 2.13 and 2.14, %(z,t) Is defined as the displacement from the mean for 
each depth over time, y is a matrix of observed covariances defined by 

<^(Zj)^(z 2 )>, E is a matrix of mode shape products defined by (^(z,) ( t ) / (z 2 ), 
and n y is noise associated with displacement estimates across depth. 

j 

= 5>, (Ofy (z) + n y (2.13) 

7=1 
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y = Ex + n y (2.14) 

In Equations 2.15 and 2.16, 8x is defined as the difference in 
displacement in the x direction (x-x a ), 8y is the difference in displacement in 

the y direction (y~y o ), S is the a-priori covariance of 8x , R is the a-priori 
covariance of the noise, and P is the covariance of Sx . 


5x=(S 1 +E T R' 1 E)' 1 E T R' 1 8y (2 15) 

P=(S 1 +E r R _1 E) _1 (2 16) 

The mode energy distribution (Figure 13) was found by taking the RMS 
value of the mode amplitudes <2, (0 . 



Figure 13. Energy (RMS( a jit ) 2 ) )of 50 modes for tidal components at DVLA and 

T1 moorings. 
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Though the wave spectrum is more complex than the mode one structure, 
most of the energy distribution is present in mode one. However, there is 
noticeable energy associated with higher modes. There may be some acoustical 
variability associated with higher modes but we will be focusing on mode 1. 
Figure 13 shows that mode one in the tidal components of the waveform 
traveling from DVLA to T1 is the largest contributor of energy in the water 
column. Therefor, we would expect to see a significant amount of acoustic 
variability from the waveform to come from the mode one tidal components. 

D. INTERNAL TIDE SOURCE REGION 

Two prominent and well-documented source regions for internal tide 
generation in the Philippine Sea are the Luzon strait in the South China Sea and 
the Sibutu Passage in the Sulu Sea (Girton et al. 2011). To develop a model for 
the region, utilizing the displacement data from DVLA and T1, a wave direction 
had to be assumed. Due to the experiments proximity to the Luzon strait, that 
location was assumed to be the source region with the propagation path heading 
due east, first passing through DVLA and then, T1. This is consistent with 
horizontal current observations during the NPAL pilot study. Taking into account 
the assumed path of propagation, the angle between DVLA and T1 was 
determined to be 72° relative to the propagation path (Figure 14). 
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Figure 14. Conceptualization of propagation path and angle between DVLA and 

T1. After https://maps.google.com/maps 

/ms?msid=213694482203726319071 ■0004e56fec1a8eab64184&msa 

=0&ll=21.371244,126.035156&spn=24.636776,27.070312 . 

E. REGIONAL MODEL 

To aid in future analysis and comparison with acoustic data obtained 
during the NPAL pilot study, a regional model was developed. The model was 
developed by combining weighted least-squares fits of the mode amplitude time 
series for the diurnal and semidiurnal frequency bands (oc /; (Q) and (3„(Q)) at 
DVLA and T1 with the previously derived 50 mode shapes ( a n (z )) extrapolated 

across 5,500 meters, the depth of the water column. The wavenumbers (k) were 
multiplied by cos(72°) to project the wavenumber along the acoustic path. In 
order to get the combined diurnal/semidiurnal (Q) amplitude components ( 
h,(x,z,t)) from the amplitudes of the two waves (C, l2 (x,t)) across the 185 km 

range (X) between DVLA and T1 for each time step (/), the amplitudes 
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themselves needed to be weighted. The wave amplitudes (£, and £ 2 ) were 
linearly weighted (w 12 (jc)) in the x-direction. At x=0km, the wave amplitude 

contributions from DVLA (£,) were weighted ( M ’i) at 100%. At x=185 km, the 


wave amplitude contributions from T1 (£ 2 ) were weighted ( w 2 ) at 100%. 50% 


weighting was applied to each wave at the midpoint. Equations 2.17 and 2.18 
were used to implement the model in Matlab. 

h "2 

wv)=XX a n (z)( a n (Q) cos(kx - Qr) + (3„ (Q) sin (kx - Q.t )) 


n=/i N=n t 

Z > (x,z,t) = £ > i(x,z,t)w l (x) + £ )2 (x,z,t)w 2 (x) 


(2.17) 

(2.18) 


Of note, the model is a data driven model that was developed for the 
Philippine Sea environment. Figures 15 and 16 were generated from the model. 
Figure 15 represents the sum of modes one and two of the tidal components and 
Figure 16 represents the sum of modes one through ten. The model has the 
ability to sum across all or a subsection of 50 modes and to separate the diurnal 
and semidiurnal components or to combine both parts in order to get the total 
internal tidal displacement between DVLA and T1. 

Figure 15 shows the simplicity of the tidal amplitudes in the ocean 
between DVLA and T1 for mode one and two. Mode one accounts for the largest 
contribution of energy, of all the modes, in the water column. Figure 16 shows 
the complexity of accounting for more modes in the wave propagation between 
DVLA and T1. As more modes are taken into account, the wave structure 
becomes more complex and the total amount of energy accounted for in the 
water column increases, represented by an increase in the magnitude and spatial 
variability of tidal amplitudes. Analysis of multiple modes in comparison with 
acoustic data is significantly more complicated than comparison with mode one. 
Fortunately, mode one accounts for the majority of the energy in the water 
column and analysis of this mode was sufficient to accomplish the objectives of 
this thesis. 
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Depth (km) Depth (km) g; Depth (km) 


Mode 1-2 Diurnal Displacement between DVLA and T1 04-Apr-2009 12:00:13 



0 20 40 60 80 100 120 140 160 180 

1e 1-2 Semidiurnal Displacement between DVLA and T1 04-Apr-2009 12:00:13 



0 20 40 60 80 100 120 140 160 180 

Mode 1-2 Tidal Displacement between DVLA and T1 04-Apr-2009 12:00:13 



0 20 40 60 80 100 120 140 160 180 


Range (km) 

Figure 15. Internal tide model for the Philippine Sea, modes 1&2. Displacement is 

in meters. 
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Mode 1-10 Diurnal Displacement between DVLA and T1 04-Apr-2009 12:00:13 
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Mode 1-10 Semidiurnal Displacement between DVLA and T1 04-Apr-2009 12:00:13 



0 20 40 60 80 100 120 140 160 180 

Mode 1-10 Tidal Displacement between DVLA and T1 04-Apr-2009 12:00:13 



0 20 40 60 80 100 120 140 160 180 


Range (km) 

Figure 16. Internal tide model for the Philippine Sea, modes 1-10. Displacement is 

in meters. 
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III. ACOUSTIC ANALYSIS 


Weston R.T. Coby conducted the acoustic analysis presented in this 
chapter. 

The data collected during the NPAL pilot study were used for the 
examination of eigenray travel time, and intensity trends at multiple depths. A 
comparison of the eigenray travel time fluctuations with the internal tide 
displacements was used to determine if there was a correlation between the 
strength of the internal tides and travel time variation. 

To conduct the acoustic analysis a prediction model was created for the 
eigenrays, and timefronts at four separate depths. The prediction model was then 
compared to the observed timefronts, and used to determine the ray path 
characteristics. Next, the Viterbi algorithm was used to track the eigenray 
timefronts. Finally, a statistical analysis of the actual acoustic data was 
conducted, to quantify the observed fluctuations. 

A. RAY THEORY AND TIMEFRONTS 

Ray theory is a high frequency approximation that describes propagation 
of sound rays through the ocean. In the deep water of the Philippine Sea, 
acoustic rays initially experience spherical spreading and later experience 
cylindrical spreading once they are vertically confined by the surface and bottom 
boundaries. As acoustic rays propagate away from a source, they are refracted, 
reflected, and scattered. This causes variations to their arrival angle, arrival time, 
and intensity of eigenrays at the receiver. Eigenrays are acoustic rays that travel 
from the transmitter and arrive at the exact location of the receiver. Ray theory 
makes the fundamental assumptions that there are no diffraction, that acoustic 
rays propagate as local plane waves, and that only small sound speed spatial 
changes take place over an acoustic wavelength. 

Over the course of the pilot study eigenrays travelled from the T1 mooring 

source at 1,095m, to a receiver at a specific depth on the DVLA. The eigenray 
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properties are described through Equations 3.1-3.6. The ray speed is described 
by the group speed Equation (3.1), which states that acoustic wave energy 
follows the ray path at the group velocity, c g . The group velocity is equal to the 
total derivative of the ray path ( dx/dt) with respect to time, and equal to the 
partial derivative of the plane dispersion relation, co (Equation 3.2) with respect 
to the vector wavenumber, k . The directional changes to the ray path are 
described by the refraction equation (Equation 3.3), which states that a changes 
to the vector wavenumber with respect to time are equal and opposite the 
changes to the dispersion relation (Colosi and Worcester 2013, unpublished 
manuscript). 


dx/ dt = 8co8k = c s 

(3.1) 

CD = (k x 2 +k y 2 +k*Y 2 c(r) 

(3.2) 

dk/dt - -8co/8x 

(3.3) 


The dispersion relation considers changes in a 3-D model (x, y, and z 
coordinates), yet the scope of this thesis only requires us to examine a 2-D ray 
slice containing the range and depth (x and z coordinates). The ray equations 
(Equations 3.4-3.6) use range as the independent variable (fixed transmitter and 
receiver) for a ray propagating in only the positive x-direction. From the CTD 
generated sound speed profiles we were able to get c(z), which we then used to 
solve the ray equations. The ray equations utilized data from both CTD profiles 
and buoys to predict the ray propagation paths from the dispersion relation, 
wavenumber (k z ), depth (z), and sound speed (c). The ray equations describe 
how the ray path changes as the ray travels from its source to the receiver. 
Because the ray equations are independent of the acoustic frequency, the ray 
paths will not change if the source frequency varies. Equation 3.4 calculates the 
ray’s refraction describing the evolution of the ray’s direction of propagation as it 
moves away from the source. Equation 3.5 calculates the ray’s change in depth 
as it propagates. Equation 3.6 calculates the change in the ray’s travel time as it 
propagates (Colosi and Worcester 2013, unpublished manuscript). 
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This thesis focuses on sound rays propagating about the sound channel 
axis (depth at which the lowest sound speed occurs). The Ray Identification 
number (ID) describes the number of turning points within a ray. Acoustic waves 
with larger launch angles have lower IDs than acoustic waves with smaller 
launch angles. In deep water, lower IDs are more easily excited, meaning that 
these rays have more energy, and that there is a greater distribution of energy 
between launch angles. As sound rays propagate they are refracted toward 
regions with the slowest sound speed. This ray refraction is observed by the 
bending of the sound rays towards the sound channel axis as they propagate 
through the water column. The rays with the steepest launch angles and lowest 
IDs are expected to arrive at the receiver first, because they travel the farthest 
from the sound channel axis (SCA) and experience higher sound speeds. 
However, rays with shallow launch angles and large IDs are expected to 
propagate closer to the SCA, experience the slowest sound speeds, and arrive 
last. 

The timefront plot (Figure 17) shows the pattern of a wavefront sweeping 
by a vertical set of receivers at a fixed range. Each wavefront branch is seen as 
the line between caustics (acoustic focusing points where the ray folds), and is 
composed of eigenrays with the same ID. The caustics farthest to the left of the 
timefront plot represent the rays with large launch angles and the fastest travel 
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times, while the caustics towards the right represent the rays with smaller launch 
angles and slower travel times. 


Timefront predictions 



Figure 17. Plot of the predicted timefront. 


The NPAL pilot study utilized a transmitter close to the sound channel axis 
to minimize acoustic interaction with boundaries, and to provide the greatest 
separation between the timefront branches (Figure 17). 

The ray trace diagram shown in Figure 17 was used to determine the 
predicted Ray IDs, and turning depths. Rays that leave the source at positive 
angles with reference to the horizontal have positive ID values, while rays with an 
initial negative angle have negative ID values. A ray with a negative ID value will 
generally arrive ahead of a ray with the same positive ID value, and nearly equal 
launch angle. 
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Eigenrays 



Figure 18. Plot of the predicted eigenrays for the 1,385m, receiver, derived from 
CTD sound speed profile. As the launch angle decreases the IDs 

increase. 


The predicted Ray IDs correspond to the number of wavefronts seen in 
the timefront plot. The ray propagating closest to the sound channel axis (Ray ID 
+9) arrives last, while the ray with the steepest launch angle (Ray ID -6) arrives 
first. 

B. ACOUSTIC RAY TRAVEL TIMES 

The pilot study consisted of 1,543 transmissions from the T1 mooring to 
the DVLA over a one-month period. The transmissions were sent out at intervals 
of 6 minutes and 3-hours. The acoustic data used for this thesis were collected 
from 24 receivers on the DVLA located at depths ranging from -800 m to 
-1385m. Four focus depths were chosen to provide insight on how the acoustic 
signals varied above the SCA (800m), along the SCA (1,000m), and between two 
receivers located within 100m of each other below the SCA (1,300m, and 
1,385m). 

An acoustic signal transmitted from a source (located at 1,095m) on the 
T1 mooring close to the SCA travelled to the DVLA mooring along several 
acoustic ray paths (Colosi et al. 2012). We used the ray equations to predict the 
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number of eigenrays, their paths, and timefronts for the four focus depths. Using 
inputs for the receiver depth, minimum and maximum launch angles, and launch 
angle step size we were able to create ray path predictions. Increasing the 
launch angle range, and decreasing the step size between launches, resulted in 
a noticeable increase in the number of predicted eigenrays. The optimum launch 
angle range was from -16° relative to the horizontal (SCA), to +16°, and the 
optimal launch angle step was .5°. Utilizing sub-optimal values decreasing the 
launch angle range, or increasing the step size resulted in fewer predicted 
eigenrays, excluding higher Ray IDs. Conversely increasing the launch angle 
range, or decreasing the step size made processing times unmanageably long. 

The number of eigenrays predicted was dependent on the depth of the 
receiver relative to the transmitter. The receiver above the sound channel (800m) 
had eight eigenrays, none of which travel along the SCA. The number of 
eigenrays increased to thirteen at 1,000m (Figure 19) with several of those rays 
travelling along the SCA. 


Eigenrays 
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Figure 19. Predicted eigenray plot for 1,000m depth receiver. Improved eigenray 
resolution revealed additional eigenrays along sound channel axis. 
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Ten eigenrays were predicted at 1,300m and 1,385m (Figure 18), losing 
the larger numbered Ray IDs that were found at the 1,000m receiver depth. 
When the number of eigenrays increases, the ability to distinguish between each 
acoustic ray decreases due to the proximity of their arrivals. This congestion of 
arrival times is most pronounced near the SCA as the four slowest eigenrays 
arrive within .005 seconds of each other (Table 2). 


Receiver Depth 

1000m 

Ray ID 

-9 

9 

9 

-9 

Arrival # 

10 

11 

12 

13 

Launch Angle (Deg) 

-4.75 

0.25 

0.75 

-1.125 

Turning Depth (m) 

652.28 

1000.23 

946.91 

910.90 

End Time (s) 

124.94 

124.94 

124.94 

124.94 


Table 2. 1,000m Receiver Ray arrival times. 


The acoustic data collected at DVLA were used to calculate SNR, which 
was then used to track the acoustic timefronts. The travel times for the tracked 


timefronts at the focus depths were then compared to the predicted arrival times 
to determine Ray IDs (Table 3). 


Receiver Depth 

800m 

Predicted ID 

-6 

-6 


6 

7 

-7 

-8 

8 

9 


Tracked Ray ID 

-6 

-6 

6 

7 

-7 


8 

9 


Receiver Depth 

1000m 

Predicted ID 

-6 | -6 6 6 

7 

-7 

-8 8 | 9 -9 | 9 9 -9 

Tracked Ray ID 


7 

-7 





Receiver Depth 

1300m 

Predicted ID 


-6 

6 

6 

7 

-7 

8 

-8 | 9 -9 | 9 


Tracked Ray ID 

-6 

6 

6 

7 

-7 

8 




Receiver Depth 

1385m 

Predicted ID 


-6 

6 

6 

7 

-7 

8 

-8 | 9 -9 | 9 


Tracked Ray ID 

-6 

6 

6 

7 

-7 

8 



Table 3. Predicted Ray IDs and tracked Ray IDs for focus depths. 
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c. 


TIMEFRONT TRACKING 


Figure 20 shows an overlay of the predicted timefronts for the receiver 
depths, and the observed timefronts for the first transmission. The timefront 
overlay plot clearly shows that for this transmission while the lower Ray IDs are 
not accurately predicted, the higher Ray ID timefronts closely conform to the 
predicted timefront. Figure 20 shows that the observed timefronts for Ray IDs 
±6 —> ±8 are easily distinguishable and relatively well spaced for the focus 
depths, while the ±9 timefronts are difficult to isolate. Furthermore, it shows that 
as depth changes the spacing between the wavefronts also varies causing some 
rays to move closer together and others to move farther apart. 


Timefront on 07-Apr-2009 06:00:00 



124 124.2 124.4 124.6 124.8 125 

Travel time (sec) 


Figure 20. Observed transmission timefront with predicted timefront overlay. 

Ray arrival plots (Figures 21-24) show the eigenray arrival times over the 
duration of the NPAL pilot study. The black dots in Figures 22-23 show the peak 
arrivals at the 1000m receiver depth with two different dB cutoffs (20dB and 
25dB). Figure 24 shows the Viterbi tracked arrivals at 1,000m. Peak arrivals and 
the Viterbi tracking method will be discussed more in depth later in the section. 
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Figures 21 through 23 show that for the receiver at 1,000m the Ray IDs 
are harder to distinguish, and Figure 24 shows that only Ray IDs ±7 can be 
tracked. What cannot be seen from the timefront plots in Figures 21 through 23 is 
that Ray IDs ±6, and ±9 each have multiple eigenrays that cannot be readily 
distinguished. 


1,300m Arrivals 
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Travel time (seconds) 

Figure 21. Time series plot of eigenray arrivals at 1,000m depth receiver. 

The ray timefronts were tracked over the duration of the experiment in 
order to conduct a statistical analysis of the eigenrays. The Matlab findpeaks 
function was used on the SNR data to isolate receptions above a given 
threshold. Through trial and error it was found that a 20dB or 25dB threshold 
provided the best ray track at each depth (Figures 22 and 23), while a 30dB 
threshold left too many gaps in the arrival information, and a 15dB threshold did 
not eliminate enough noise. A 25db threshold limit was used in this thesis 
because it provided the cleanest eigenray tracks for further analysis. 
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Travel time (seconds) 

Figure 23. Time series plot of eigenray arrivals at 1,000m depth receiver with 

25dB threshold SNR peaks 

The Viterbi program is an automated tracking algorithm that was used on 
the 25db data to find the best track fit for each of the eigenrays (Dzieciuch and 
Worcester 2013, unpublished manuscript). The Viterbi program was designed to 
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Figure 22. Time series plot of eigenray arrivals at 1,000m depth receiver with 

20dB threshold SNR peaks 
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track the timefront for a single acoustic ray. The Viterbi program accuracy was 
dependent on the ability to manually set the appropriate margins of error for each 
individual timefront. The program was capable of accurately tracking well-spaced 
timefronts with a high SNR, yet was not able to track closely spaced timefronts. 
The close proximity of the arrivals at the 1,000m receiver allowed only Ray IDs ±7 
to be tracked by the Viterbi program (Figure 24). The Viterbi program works by 
creating a window, encompassing multiple peaks (Figure 23) for each sampling 
period, that is likely to contain the correct track. The Viterbi program then uses a 
Kalman filter to find the next most likely ray track position (Dzieciuch and 
Worcester 2013, unpublished manuscript), vice tracking only the highest SNR 
positions, or creating the straightest track from the first to the last reception. To find 
the most likely path, the Viterbi program utilized multiple physical parameter inputs, 
and limits, which required adjustments for each eigenray to create the most 
accurate track. 



Figure 24. Time series plot of eigenray arrivals at 1,000m depth receiver with 

Viterbi track overlay 
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The ability to track each individual eigenray depended on correcting the 
allowable standard deviation (STD) for the next peak arrival time, and the STD 
next peak observation parameters in the Viterbi program. The less critical 
parameters used to fine-tune the eigenray track were the probability that the 
track continued/changed, and the probability of the track entering the designated 
gate (window). The windows needed to be sufficiently small to prevent the target 
eigenray’s track from crossing over with the adjacent ray tracks. Decreasing the 
window size too much was problematic, as the ray tracks experienced large 
changes to the eigenray arrival times over a short period of time. These 
fluctuations resulted in the Viterbi program either losing the eigenray track, or 
crossing over between tracks. As noted earlier, the ray predictions for 1,000m, 
1,300m, and 1,385m contain more than eight eigenrays. However, due to 
insufficient spacing between the highest Ray IDs of the 1,300m, and 1,385m 
eigenrays, the Viterbi program was only able to track the first six eigenrays. 
Similarly, due to the proximity of the adjacent eigenrays at 1,000m, only the 
middle two eigenrays could be tracked using the Viterbi algorithm. 

The fluctuations in the ray time series (Figure 24) can be quantified by the 
travel time variation (Equation 3.7) and scintillation. The change in travel time (St) 
due to sound speed fluctuations (8c) over the ray path (s) is shown in equation 
3.7. The effect on intensity due to sound speed fluctuations is much more 
complicated and cannot be easily explained. 



D. ACOUSTIC TIMEFRONT VARIABILITY 

The ability and inability to track eigenrays over the duration of the NPAL 
pilot study led to the discovery of certain significant features. The eigenray arrival 
times varied significantly over the duration of the experiment, yet similar trends 
were observed for the tracked rays across the focus depths. 
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Over certain reception intervals the focus depths experienced similar 
drastic changes to the arrival times (Figure 25). The first 14 days of the 
experiment (07Apr 0000Z -21Apr 0300Z) are significant because of the highly 
variable SNR peaks, and eigenray arrivals, which made accurate tracking the 
individual eigenrays quite difficult. From 1145Z-1630Z on 30 April there was a 
significant increase in the eigenray travel times following a period of minimal (if 
any) SNR data, which prevented the Viterbi program from properly tracking the 
eigenray arrivals. 



Figure 25. Observed SNR trends for time series plot of eigenray arrivals at 

1,300m receiver depth (20db threshold) 


At each of the focus depths Ray ID +9 had a significantly higher SNR value 
than the other rays. However, the Ray ID associated with the lowest SNR value varied 
with depth as shown in tables 4-7. The SNR values individual rays is significant 
because it shows that the the proximity of the eigenray paths is the determinate factor 
in whether or not the Viterbi algorithm can track the ray. Figure 26 shows the inability of 
the Viterbi program to track the Ray ID -8 arrival at the 800m receiver due to its initial 
proximity to Ray ID +8. 


39 










Figure 26. Time series plot of eigenray arrivals at 800m depth receiver with Viterbi 

track overlay. 

E. TRAVEL TIME AND INTENSITY STATISTICS 

A statistical analysis was conducted on each of the tracked timefronts to 
better understand how the ray timefronts varied throughout the experiment. The 
mean arrival time, arrival time variance, arrival time STD (wander), and the errors 
were calculated for each eigenray at the focus depths. Additional calculations 
were done to determine the peak broadband Intensity (I), Scintillation Index (SI) 
(Equation 3.8), intensity STD, and SI error. While the intensity calculations are 
too difficult to briefly explain, the SI is a measure of the intensity fluctuations. 

(/ 2 ) 

SI = \pf~ 1 (3 ' 8) 

Tables 4-7 contain the data from the SI, and the statistical analysis for SNR 
and time. Tables 4-7 show that the eigenrays with the greatest launch angles arrive 
first yet have much lower SNR values than the eigenrays with the shallowest launch 
angles. The SNR variation for a fully saturated environment is approximately 5.6dB. 
From the tables it appears that only the rays propagating closest to the sound 
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channel axis experience SNR variation that is consistent with a saturated 
environment. As would be expected the scintillation index also appears to generally 
increase with decreasing launch angles. From the tables the phase variability 
appears to generally increase for all eigenrays as they approach the SCA and 
decrease as they move away from the sound channel axis. 


Receiver Depth 

800m 

Ray ID 

-6 

-6 

6 

7 

-7 

8 

9 

Arrival # 

1 

2 

3 

4 

5 

7 

8 

Launch Angle (Deg) 

-14.72 

-15.14 

14.59 

13.38 

-10.25 

9.03 

8.09 

Turning Depth (m) 

42.51 

0.02 

58.43 

171.36 

412.13 

468.96 

508.73 

Arrival Mean (s) 
Arrival STD (s) 
Track Variance (s) 

123.99 

124.10 

124.25 

124.64 

124.71 

124.82 

124.90 

1.18E-02 

1.65E-02 

2.1 IE-02 

1.82E-02 

1.50E-02 

1.84E-02 

1.83E-02 

1.38E-04 

2.71E-04 

4.46E-04 

3.3 IE-04 

2.25E-04 

3.37E-04 

3.34E-04 

SNR Mean (dB) 
SNR STD (dB) 

32.38 

34.09 

30.15 

31.60 

31.61 

32.91 

39.51 

3.30 

4.50 

3.38 

3.65 

3.85 

4.01 

6.52 

Scintillation Index 

2.79 

4.69 

5.95 

4.66 

4.24 

3.16 

3.97 

Phase Variability (rad/s) 

20.31 

28.44 

36.48 

31.44 

25.91 

31.74 

31.59 


Table 4. 800m Receiver tracked timefront statistics and SI. SNR statistics 

account for internal wave and mesoscale variability. 


Receiver Depth 

1000m 

Arrival # 

3 

4 

Ray ID 

7 

-7 

Launch Angle (Deg) 

13.28 

-10.38 

Turning Depth (m) 

180.85 

405.73 

Arrival Mean (s) 

124.28 

124.62 

Arrival STD (s) 

2.24E-02 

1.80E-02 

Track Variance (s) 

5.00E-04 

3.25E-04 

SNR Mean (dB) 

30.84 

30.74 

SNR STD (dB) 

3.26 

3.64 

Scintillation Index 

0.494 

0.578 

Phase Variability (rad/s) 

38.65 

31.13 


Table 5. 1,000m Receiver tracked timefront statistics and SI. SNR statistics 

account for internal wave and mesoscale variability. 
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Receiver Depth 

1300m 

Arrival # 

3 

1 

2 

4 

5 

6 

Ray ID 

-6 

6 

6 

7 

-7 

8 

Launch Angle (Deg) 

-14.59 

14.72 

15.14 

13.13 

-10.55 

9.41 

Turning Depth (m) 

58.43 

42.51 

0.02 

194.90 

396.68 

452.30 

Arrival Mean (s) 

124.33 

124.00 

124.08 

124.58 

124.71 

124.76 

Arrival STD (s) 

2.19E-02 

1.54E-02 

1.34E-02 

1.72E-02 

1.92E-02 

1.67E-02 

Track Variance (s) 

4.78E-04 

2.38E-04 

1.78E-04 

2.97E-04 

3.67E-04 

2.79E-04 

SNR Mean (dB) 

31.18 

31.68 

33.12 

30.77 

31.87 

37.51 

SNR STD (dB) 

3.32 

3.14 

5.02 

3 5H 

3 53 

6.05 

Scintillation Index 

0.478 

0.421 

2.032 

0.526 

0.596 

1.621 

Phase Variability (rad/s) 

37.78 

26.65 

23.07 

29.78 

33.11 

28.85 


Table 6. 1,300m Receiver tracked timefront statistics and SI. SNR statistics 

account for internal wave and mesoscale variability. 


Receiver Depth 

1385m 

Arrival # 

3 

1 

2 

4 

5 

6 

Ray ID 

-6 

6 

6 

7 

-7 

8 

Launch Angle (Deg) 

-14.56 

14.74 

15.12 

13.09 

-10.59 

9.47 

Turning Depth (m) 

61.06 

39.99 

0.00 

198.64 

394.16 

449.47 

Arrival Mean (s) 

124.34 

123.99 

124.09 

124.57 

124.70 

124.77 

Arrival STD (s) 

2.22E-02 

1.49E-02 

1.31E-02 

1.70E-02 

1.97E-02 

1.62E-02 

Track Variance (s) 

4.93E-04 

2.22E-04 

1.72E-04 

2.90E-04 

3.90E-04 

2.63E-04 

SNR Mean (dB) 

31.35 

31.38 

35.55 

30.67 

31.92 

37.75 

SNR STD (dB) 

3.49 

3.34 

5.06 

3.53 

3.83 

5.61 

Scintillation Index 

0.559 

0.433 

0.861 

0.603 

0.679 

1.029 

Phase Variability (rad/s) 

38.38 

25.77 

22.67 

29.43 

34.11 

28.00 


Table 7. 1,385m Receiver tracked timefront statistics and SI. SNR statistics 

account for internal wave and mesoscale variability. 


Figure 27 is a plot of the arrival wander for each of the eigenrays at the 
focus depths. The plot appears to show a general increase in arrival wander as 
the travel time increases. 
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Figure 27. Plot of eigenray wander across the focus depths. 


Figure 28 is a plot of the eigenray intensity variance across the focus 
depths. This plot shows that there is a general increase in the intensity variance 
as you increase travel time. 
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Figure 28. Plot of eigenray intensity variance across the focus depths. 
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IV. OCEAN DATA AND ACOUSTIC DATA COMPARISON 


Jacob A. Fischer, and Weston R.T. Coby jointly conducted the acoustic 
and oceanographic data comparison presented in this chapter. 

A. TIDAL FIT TO ARRIVAL TIMES 

Having determined that the internal tidal frequencies were the dominant 
frequencies, that mode 1 had the most energy within the internal wave traveling 
between DVLA and T1, and having successfully tracked individual rays through 
data analysis and using the Viterbi algorithm, the next step was to determine how 
much of the acoustic arrival time variability could be attributed to the internal 
tides. First, we needed to verify that there was a strong diurnal and semidiurnal 
signal in the acoustic arrival time fluctuations. This was accomplished by fitting 
the previously mentioned diurnal and semidiurnal frequencies to the Viterbi 
tracked ray arrival times by using the same weighted least-squares method 
discussed in Section ll-B of this thesis. Data error in performing the tidal fit was 
assumed to be 5 milliseconds (1/sampling rate). As seen in Figure 29, there is a 
distinct diurnal and semidiurnal frequency pattern in the arrival time fluctuations 
for individual ray arrivals at each sensor depth. 
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Tidal Fit for Ray ID +6, 800m 



Ray ID +6, 1300m 



Ray ID +6, 1385m 



Figure 29. Internal tidal fit to the +6 Ray ID for sensor depth 800m, 1300m, and 

1385m. 

The sampling frequency changed during the monthlong study. Areas with 
a higher concentration of black dots in Figure 29, corresponds with times when 
the sampling frequency was 10/hour. Areas with a lesser concentration of black 
dots, corresponds with times when the sampling frequency was 8/day. There is a 
better fit during periods where the sampling frequency was 10/hour, as compared 
to other periods where the sampling frequency was 8/day, because there were 
more observations during those times. This brings up speculation that the fit may 
be unreliable in areas where sampling frequency was 8/day. 

To evaluate the fit and determine the percentage of acoustic arrival time 
fluctuations that could be attributed to the tidal frequencies, a comparison of 
variability was conducted using Equation 4.1 for each of the tracked rays across 
the monthlong time period. 
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(4.1) 


'T'T? _ ^ 

1 r % - ~Zi 

®Arr(t) 

The percent acoustic arrival time fluctuation variability attributed to the 
internal tides ( TF % ) for the entire month is between 20% and 70%, depending on 

receiver depth and Ray ID (Figure 30). 
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Figure 30. Acoustic arrival time variability attributed to internal tides for all tracked 

rays at each sensor depth. 


To test the fit during periods of infrequent sampling, TF % was calculated 

where sampling frequency was low, 8/day. We determined that the fit was slightly 
under predicting some of the arrivals by as much as 30% with a range of 0%- 
40%. Testing the fit during periods when the frequency sampling was high, 
10/hour, 77^was between 30-70%. These are expected results due to the least- 
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squares fit having more data/unit time during the higher frequency sampling 
periods. Taking the analysis one step further, we fit only the dominant diurnal 
(K-i) and dominant semidiurnal (M 2 ) frequencies. This showed a TF % drop of 

10%. This is an expected result due to the acoustic fluctuations being effected by 
multiple tidal frequencies; not just the two dominant frequencies. 

For further comparison, the tidal fit of the acoustic arrival time fluctuations 
for Ray ID +6 at receiver depth 800m was overlaid on the mode 1 tidal fits of the 
oceanographic data for both DVLA and T1 (Figure 31). As was previously 
determined, there is a strong semidiurnal signal at T1 and a strong diurnal signal 
at DVLA. Two intriguing time windows to look at on this figure are between April 
12-22 on the top panel and April 22-May 7 on the bottom panel. It is clear that the 
acoustic fluctuations match well to the mode 1 signal at DVLA between April 12- 
22 and to the mode 1 signal at T1 between April 22-May 7. From this, we could 
glean that the acoustic variability of Ray ID +6 at receiver depth 800m was 
affected more by the internal tidal signal near DVLA during the April 12-22 
timeframe than by the signal near T1. However, during the April 22-May 7 
timeframe, the mode 1 component of the internal tide near T1 had a greater 
impact on acoustic variability. Further analysis would be required to verify this 
hypothesis. 
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Mode 1 Internal Tide Fit for DVLA and Arrival Time Fluctuations (Ray ID +6, 800m rcvr depth) 
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Mode 1 Internal Tide Fit forTI and Arrival Time Fluctuations (Ray ID +6, 800m rcvr depth) 



Figure 31. Tidal fit overlays for Oceanographic (black line) and Acoustic (red line) 
data sets. (Oceanographic data - Internal tide fit to Mode 1 amplitudes 
at DVLA and T1. Acoustic data - Internal tide fit to arrival time 

fluctuations) 

B. EFFECTS OF SPRING AND NEAP CYCLES 

Upon determining that acoustic arrival time fluctuations were influenced by 
internal tides, the next step was to determine the times during the experiment 
when tidal amplitudes were strongest and weakest. This period of time is known 
as the tidal range, or more commonly, the spring/neap cycle. Tidal range is 
determined by the phase of the two dominant semidiurnal frequencies, M2 and 
S2. When M2 and S2 are in phase, you will have the strongest tidal amplitude 
and when out of phase, the weakest. 
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The spring and neap cycles were found by utilizing the previously fit 
internal tide displacement data. 50 mode amplitudes were fit to the data across 
time using the method described in Section ll-C of this thesis. However, instead 
of fitting 15 frequencies to the mode amplitudes, only two frequencies, the M2 
and S2 were fit. From the mode one time series, times corresponding to 
maximum and minimum amplitudes were identified. The inflection point between 
the maximum and minimum displacement amplitudes was used as the dividing 
time between spring and neap cycles. The spring/neap cycle where we had the 
most complete overlap of acoustic and oceanographic data, from 13 April, 2009 
to 27 April, 2009, was chosen for analysis. Figures 32-33 show arrival time 
fluctuations overlaid on the mode displacement amplitudes. You can see from 
Figures 32-33 that there is some acoustic arrival time variability associated with 
the spring and neap cycles. 
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DVLA Mode 1 Internal Tide Fit and Ray ID +6 fluctuations for 800m receiver depth 
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DVLA Mode 2 Internal Tide Fit and Ray ID +6 fluctuations for 800m receiver depth 
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Figure 32. Mode 1 & 2 time series of amplitudes fit to the internal tides at DVLA 
shown with a red line. Arrival time fluctuations for Ray ID +6 at the 
800m receiver depth are marked with black asterisks. 
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T1 Mode 1 Internal Tide Fit and Ray ID +6 fluctuations for 800m receiver depth 
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T1 Mode 2 Internal Tide Fit and Ray ID +6 fluctuations for 800m receiver depth 



Figure 33. Mode 1 & 2 time series of amplitudes fit to the internal tides at T1 
shown with a red line. Arrival time fluctuations for Ray ID +6 at the 
800m receiver depth are marked with black asterisks. 

It can be seen by looking at Figure 32 that arrival time fluctuations for Ray 
ID +6 at the 800m sensor depth matches well with the mode 1 amplitude time 
series for DVLA. This highlights the connection between arrival time fluctuations 
and internal tide mode structure and also, gives credence to the foregone 
conclusion that the acoustic travel time fluctuations for this ray were affected by 
the semidiurnal component of the internal tide between DVLA and T1. However, 
the degree to which acoustic arrival time fluctuations are effected by the modal 
structure is a complicated puzzle, as can be seen in Figure 33, where there 
appears to be minimal connection with the mode 2 structure at T1. Comparison 
of modal structure to arrival time fluctuations is beyond the scope of this thesis 
but analysis of this topic may prove beneficial in future work. 
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To further analyze the Viterbi tracked ray arrivals during the spring and 
neap times, we needed to verify that the difference between the mean arrival 
times for each part of the cycle was statistically significant. Since the sample 
sizes and variances were different, a two-sided T-test was performed on each 
data set. An alpha of 0.01 was used. It was determined the difference between 
the means of the spring and neap time periods for all of the tracked rays, proved 
to be, with 99% confidence, statistically significant. The next step was to 
calculate the VAR of the arrival time fluctuations for each tracked ray during the 
spring and neap cycle (Figure 34). 


x 10 4 Tracked rays spring/neap cycle variance, Sensor Depth 800m 



Ray ID 


Figure 34. Tidal range variance for all tracked rays at each sensor depth. 

The variances in Figure 34 show that every ray we tracked was affected 
by the spring/neap cycle. This is significant because it shows that arrival time 
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fluctuations are directly tied to internal tidal amplitude fluctuations. For larger tidal 
amplitudes, you would expect to see more variability in arrival time fluctuations. 

C. FUTURE WORK 

Compared to the follow-on 2010-2011 NPAL experiment the data 
collected during the NPAL pilot study were limited in scope, and data resolution. 
The NPAL 2010-2011 study was conducted over the course of a year, and had 
greater time and depth resolution than the 2009 pilot study. 

The data from the 2010-2011 experiment should be used to track the 
timefronts across all depths, and study the seasonal low frequency variability due 
to eddies and other mesoscale disturbances. The data should also be used to 
resolve more tides. A comparison of tracked peaks with the internal tide model 
for multiple spring/neap cycles can then be used to provide a more 
comprehensive understanding of how evolving tides affect acoustic travel times 
(Worcester, 2008). 

A more accurate modal analysis can be conducted on the 2010-2011 data 
with the increased number of sensors distributed throughout the water column, 
requiring less data extrapolation for depth than the 2009 data. The modal 
analysis conducted within this thesis showed that modes higher than mode 1 
possessed a significant amount of energy, contrary to previous conclusions. The 
data from the NPAL 2010-2011 study should be used to compare acoustical data 
for higher internal tide modes in this complex region, to better understand the 
effect of the higher modes on acoustic propagation. 

Finally, data from the DVLA and the other moorings should be used to 
determine if there is a difference in the travel time variation depending on the 
bearing angle from the source to the receiver. This will require a comparison of 
the phase changes for acoustic signals propagating in reciprocal directions. 
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